package
{
	import flash.utils.ByteArray;
	
	/* TODO: transform unknown residues into '*' */
	
	public class Compute
	{
		/* http://en.wikipedia.org/wiki/Needleman-Wunsch_algorithm */
		/* http://snippets.dzone.com/posts/show/2199 */
		
		public static var Alphabet:Array;
		public static var SimMatrix:Array;
		
		/* Parses a similarity matrix, such as BLOSUM62. */
		public static function LoadSimMatrix(buffer:ByteArray):void {
			var MatrixLines:Array = buffer.toString().split("\n");
			var MatrixRaw:Array = new Array();
			for each (var str:String in MatrixLines) {
				var com:int = str.indexOf("#"); /* comment in this line? */
				if (com != -1) str = str.substring(0, com); /* remove comments */
				if (!str.match(/\S/)) continue; /* skip empty lines */
				MatrixRaw.push(Utils.Tokenize(str));
			}
			SimMatrix = new Array();
			Alphabet = MatrixRaw.shift(); /* First row is the alphabet */
			for each (var sm:String in Alphabet) SimMatrix[sm] = new Array();
			for each (var sym:String in Alphabet) {
				var row:Array = MatrixRaw.shift();
				if (row.shift() != sym) throw new Error("Bad Matrix!");
				for each (var rsym:String in Alphabet) SimMatrix[rsym][sym] = parseInt(row.shift()); /* Allow floats? */
			}
			debug("Similarity Matrix Loaded!");
		} /* TODO: proper error checking? */
		
		
		/* Sequence Alignment of the given sequences, using the current similarity matrix. */
		public static function NeedlemanWunsch(seq:Array, ref:Array, gap:int):Object { /* todo: return result */
			var rows:uint = seq.length + 1;
			var cols:uint = ref.length + 1;
			var a:Array = Utils.MDArray(cols, rows);
			for (var i:uint = 0; i < rows; i++) a[i][0] = gap * i;
			for (var j:uint = 0; j < cols; j++) a[0][j] = gap * j;
			for (i = 1; i < rows; i++) {
				for (j = 1; j < cols; j++) {
					var choice1:int = a[i-1][j-1] + SimMatrix[ref[j-1]][seq[i-1]];
					var choice2:int = a[i-1][j] + gap;
					var choice3:int = a[i][j-1] + gap;
					a[i][j] = Math.max(choice1, choice2, choice3);
				}
			}
			//Utils.DumpMatrix(a);
			var rf:Array = new Array();
			var sq:Array = new Array();
			i = seq.length; // rows
			j = ref.length; // cols
			while (i > 0 && j > 0) {
				var score:int = a[i][j];
				var score_diag:int = a[i-1][j-1];
				var score_up:int = a[i-1][j];
				var score_left:int = a[i][j-1];
				var sm:int = SimMatrix[ref[j-1]][seq[i-1]];
				//debug("i="+i+" j="+j+" score="+score+" diag="+score_diag+" up="+score_up+" left="+score_left+" sm="+sm);
				switch (score) { /* Counter-Clockwise */
					case (score_up + gap):
						//debug("MOVE UP");
						rf.unshift("-");
						sq.unshift(seq[i-1]);
						i--;
						break;
					case (score_diag + sm):
						//debug("MOVE DIAG");
						rf.unshift(ref[j-1]);
						sq.unshift(seq[i-1]);
						i--;
						j--;
						break;
					case (score_left + gap):
						//debug("MOVE LEFT");
						rf.unshift(ref[j-1]);
						sq.unshift("-");
						j--;
						break;
				}
				//debug(" ");
			}
			while (j > 0) {
				rf.unshift(ref[j-1]);
				sq.unshift("-");
				j--;
			}
			while (i > 0) {
				rf.unshift("-");
				sq.unshift(seq[i-1]);
				i--;
			}
			return {A:sq, B:rf}; /* Aligned sequences */
		}
		
		public static function AlignTester():void {
			debug("Test Align Started.");
			var res:Object = NeedlemanWunsch(Tests.TestSeqA, Tests.TestSeqB, -10);
			debug("A:"+res.A.toString());
			debug("B:"+res.B.toString());
			debug("Length of Output A = "+res.A.length);
			debug("Length of Output B = "+res.B.length);
			debug("Test Align End.");
		}
				
		private static function debug(msg:String):void {
        	np.debug.print(msg);
        }
	}
}